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Abstract 


A two-phase, mass-transport model for liquid-feed direct methanol fuel cells (DMFCs) is developed by taking into account the effect of non- 
equilibrium evaporation and condensation of methanol and water. The comparison between the present model and other models indicates that 
the present model yields more reasonable predictions of cell performance. Particularly, it is shown that the models that invoke a thermodynamic- 
equilibrium assumption between phases will overestimate mass-transport rates of methanol and water, thereby resulting in an inaccurate prediction 
of cell performance. The parametric study using the present model reveals that the gas coverage at the flow channel-—diffusion-layer interface is 
directly related to the gas-void fraction inside the anode porous region; increasing the gas-void fraction will increase the mass-transfer resistance 
of methanol and thus lower cell performance. The effects of the geometric dimensions of the cell structure, such as channel width and rib width, 
on cell performance are also investigated with the model developed in this work. 


© 2007 Elsevier B.V. All rights reserved. 
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1. Introduction 


The direct methanol fuel cell (DMFC) has been regarded 
as an attractive power source for portable and microelectronic 
devices because of its high efficiency, high specific energy, low 
emissions, and simple structure [1]. Nevertheless, widespread 
commercialization of DMFC technology is still hindered by 
several technological challenges, such as the low rate of the 
methanol oxidation reaction (MOR) on the anode, and mixed 
overpotential caused by methanol crossover to the cathode. In 
addition, coupled electrochemical reactions and transport of dif- 
ferent species, e.g., reactants, products, electrons, protons and 
heat, make it difficult to achieve optimum design and opera- 
tion of the DMFC. In an effort to resolve these challenging 
problems, mathematical modelling of the DMFC system plays 
an important role, as it can provide a powerful and economi- 
cal tool to analyze the complex transport processes, which are 
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rather difficult to be studied experimentally. Most of the previ- 
ous fuel cell models have been based on the single-phase flow 
approach [2-17], which is apparently unrealistic because the 
coexisting liquid and gas flow behaviour has a significant influ- 
ence on the mass-transport processes, as revealed by in situ flow 
visualization studies on DMFCs [18-20]. 

Although more realistic, the development of a two-phase, 
mass-transport model for the DMFC is rather challenging 
because of the complexity of two-phase flows. Murgia et al. 
[21] presented a one-dimensional model based on phenomeno- 
logical transport equations for the DMFC. In order to consider 
two-phase flow interaction in the diffusion layer, they introduced 
a Gaussian function as an approximate representation of influ- 
ence of the capillary pressure on the effective gas porosity. Wang 
and Wang [22] modelled the DMFC using the mixture multi- 
phase flow model, in which the species in the liquid and gas 
phase are assumed to be in thermodynamic equilibrium. Birg- 
ersson et al. [23] studied numerically the anode mass transport 
in two-phase flow in a two-dimensional domain using the math- 
ematical model developed by Wang and Wang [22]. Divisek 
et al. [24] presented a two-dimensional, two-phase model for 
the DMFC with the computational domain restricted to the 


W.W. Yang, T.S. Zhao / Journal of Power Sources 174 (2007) 136-147 137 


Nomenclature 


Reontact 
S 
T 


interfacial specific area between liquid and gas 
phase (m? m7?) 

specific area (m? m~?) 

molar concentration (mol m7?) 

diffusivity (m? s7!) 

Faraday constant (96,478 C mol`!) 

current vector (A m7?) 

current density (A m`?) 

cell current density (A m7?) 

parasitic current resulting from methanol 
crossover (A m7) 

exchange current density (A m~?) 

anode current density (A m~?) 

cathode current density (A m7?) 

condensation rate (mol atm—! s7! m7?) 
evaporation rate (1 atm—! s7!) 

Henry’s law constant (Pa) 

relative permeability 

permeability of porous material (m?) 

source term in mass conservation equation 
(kg m~? s7!) 

molecular weight (kg mol~!) 

electro-osmotic drag coefficient 

mol flux (mol m2 s7!) 

capillary pressure (Pa) 

gas-phase pressure (Pa) 

liquid-phase pressure (Pa) 

gas constant (J mol`! K7!) 

source term in species conservation equation 
(mol m~3 s7!) 

interfacial species transfer rate (mol m~? s7!) 
Ohmic contact resistance (Q m7?) 

liquid saturation 

temperature (K) 

thermodynamic-equilibrium voltage (V) 

cell voltage (V) 

coordinate, m, or mole fraction in liquid solution 
(mol mol7!) 

coordinate, m, or mole fraction in gas mixture 
(mol mol`!) 


Greek letters 


gas-void fraction (1 — s) 

anode transfer coefficient at anode 

cathode transfer coefficient at cathode 
reaction order of ORR 

thickness of porous layer (m) 

porosity of porous medium 

overpotential (V) 

contact angle (°) 

ionic conductivity of membrane (Q7! m7!) 
viscosity (kgm7! s7!) 


p density (kg m~?) 

o interfacial tension (N m7!) 
Subscripts 

A anode 

Acl anode catalyst layer 

Adl anode diffusion layer 

Afc anode flow channel 

C cathode 


Ccl cathode catalyst layer 
Cdl cathode diffusion layer 


Cfc cathode flow channel 
G gas phase 

L liquid phase 

Mem membrane 

M methanol 

MV methanol vapour 

W water 


WV water vapour 


Superscripts 

eff effective value 

in inlet condition 

ref reference value 
sat saturated value 


membrane electrode assembly only. In their work, the presence 
of hydrophilic and hydrophobic pores was taken into account 
by introducing an experimental-data-fitted correction equation 
between capillary pressure and saturation, which was differ- 
ent from the Leverette function used in Ref. [22]. Although 
methanol permeation was considered, the mixed-potential effect 
was not addressed. Rice and Faghri [25] proposed a tran- 
sient, multiphase model for a passive fuel-fed DMFC. The 
passive delivery system consisted of a porous material to deliver 
methanol to the fuel cell while controlling the concentration 
of methanol at the anode. This procedure allowed the rate of 
methanol crossover to be reduced. Noteworthy is that the evap- 
oration/condensation of methanol and water was formulated in 
a manner to capture non-equilibrium effects between phases. 
This approach differed from the thermodynamic-equilibrium 
assumption employed elsewhere [22,23]. Recently, Liu and 
Wang [26] presented a two-phase model based on their ear- 
lier version [22], in which the catalyst layer was treated as a 
component, rather than an interface without thickness. In their 
DMFC model, however, the liquid and gas phase were assumed 
to be in thermodynamic equilibrium. Most recently, Yang and 
Zhao [27] presented a two-phase mass-transport model that was 
developed based on multiphase flow theory in a porous medium 
and hence there was no need to define the mixture pressure of 
the liquid and gas phases as in the multiphase mixture model. 
Another feature of the model was that the assumption of con- 
stant gas-phase pressure is eliminated. On the other hand, the 
model did not take into account the effect of methanol vapour 
transport. 
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The objective of this work is to extend the isothermal, 
two-phase, mass-transport model for DMFCs that has been 
developed earlier [27] by taking into account the effect of 
non-equilibrium evaporation and condensation between the 
phases. This new model enables removal of the assumption 
of a thermodynamic-equilibrium condition between the phases 
employed in previous models [22,23,26]. The model also elim- 
inates the assumption of constant liquid pressure at the anode 
or constant gas-phase pressure at the cathode, as was adopted 
previously [28-30]. The mathematical model is presented and 
compared with other models. Finally, numerical results obtained 
from parametric studies by the model are presented and dis- 
cussed. 


2. Mathematical model 


The two-dimensional domain, as sketched in Fig. 1, repre- 
sents the typical geometry of a membrane electrode assembly 
(MEA) that consists of anode diffusion layer (ADL), an anode 
catalyst layer (ACL), a membrane (MEM), a cathode catalyst 
layer (CCL) and a cathode diffusion layer (CDL). The MEA is 
sandwiched between two parallel flow-field plates. Since both 
channel width and rib width in the parallel flow-field are sym- 
metrical with respect to their middle point, only a half-rib width 
and a half-channel width need to be considered to save the com- 
puting time. As such, the computation domain is the region that 
is enclosed by the dashed lines. 


2.1. Governing equations 


We now present the governing equations for the two- 
dimensional, steady-state, two-phase mass transport in the 
DMFC anode and cathode porous regions. The model is devel- 
oped based on our previous work [27], in which transport of 
methanol vapour and water vapour was not taken account. To 
consider the transport of these additional species with phase 
change, the effect of non-equilibrium evaporation and con- 


Channel 


Fig. 1. Schematic of model domain. 


densation between the phases needs to be taken into account. 
This opposes the assumption of a thermodynamic equilibrium 
between the phases that was employed in previous models 
[22,23,26]. As such, in the anode and cathode porous region, a 
total of nine variables are involved. These includes the concen- 
tration of liquid methanol (Cy), the concentration of methanol 
vapour (Cmy,g), liquid water saturation in the anode (s1,a), liquid 
pressure in the anode (pi a), the concentration of water vapour 
in the anode (Cwvy,g,a), the concentration of gas oxygen in the 
cathode (Co,,g), liquid water saturation in the cathode (s),c), 
gas-phase pressure in the cathode (pg), and the concentration 
of water vapour in the cathode (Cwy,g,c). Note that methanol 
vapour concentration, Cyy,g and water vapour concentration, 
Cwv,g,a are the two added variables compared with the pre- 
vious model [27]. The governing equations that describe the 
mass conservation of each species in different phases are given 
below. 


2.1.1. Anode porous region 


Kk eff : 
CM: V- f= a Vpi | Cui — DM VCmia} = Ruyia (1) 


Kk 
CMV.g : Ve I(- le pea) CMV.,g,a = Dil gVCMV 09 


o 
5 


= RMV, g.a (2) 


Kk, d 
Sha: V. (- 2 ) ( Pe VSl,a + Ves) = Mega (3) 


Hg/Pg ds 
N. 
V Pea 
. Kky : 
Pl,a : V. — V Pla = Mj] a (4) 
Mi/ PI 


Kk, 
Cw. :V. (- A pea) CwV,g,a = Diy gVEwv.e| 


o 
(= 


= Rwy., g,a (5) 


2.1.2. Cathode porous region 


Kk 
Cosg: V> I(- le Vps) Cos,g,c — D8 VCore 
8 


= Ro,, gc (6) 


Kky d . 
Sie Ve ( r ) ( $ c Ysie + pse) =m (7) 
a, 


Kkro 
Pg,c : V. (- = ) Ype = Mec (8) 
Mg/Pg 
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_ Kkrg _ neff 
Cwv.gc: V- r V Pg,c | Cwv,g,c — Dwy,g VCwv,g,c 
g 


= Rwy ec (9) 


The interfacial transfer rates between liquid and vapour of 
methanol and water in the anode porous region are, respectively, 
by 


Ru = Aighigs( — s PM — Py (10) 
RT 
and 
ESPI P ' 
kez (wv Pe — Pw) Ywv Pe < pry 
4 1,0 
Rw (11) 
e(l — s)ywv sat sat 
ke —__ (wv Pg — Pwy) YwvPs > Pwy 


RT 


where pi, and pý, denote the saturation pressure of methanol 


vapour and the saturation pressure of water vapour, respectively. 


2.1.3. Membrane 

Unlike in the anode porous region, only dissolved methanol 
and liquid water need to be considered, in the polymer elec- 
trolyte membrane as this is usually regarded as a gas insulator. 
Transport of methanol through membranes generally depends 
on molecular diffusion, electro-osmotic drag, and convection. 
Accordingly, the flux of methanol crossover (Nm) is given by 


I K Apic: 
Nm = -DË . VCm + nam ( Eb :) Cu (12) 


(Nw) is given by 
I pı K Apic-a 
F  Mp,o Hi mem 


Nw = na, mo (13) 

To make the above governing equations closed, some con- 
stitutive correlations and definitions are needed. These include 
capillary pressure, relative permeability for both gas and liq- 
uid phases, effective diffusion coefficients for each species, 
source terms and interfacial transfer rates of methanol and 
water between the phases. All these correlations and associated 
nomenclatures are listed in Table 1. 


2.2. Boundary and interfacial conditions 


As indicated in Fig. 1, the computational domain is enclosed 
by six boundaries. The conditions at each boundary are described 
below: 


e Boundary 1: This boundary represents the inlet of reactant 
supply at the anode, at which the concentration of liquid 
methanol, liquid-phase pressure, liquid saturation and the 
concentrations of methanol vapour and water vapour are all 
specified to be inlet conditions: 

CM. = Cc 


Cuv = Cu: Cwv = Cy; 


p=p s=! (14) 

e Boundary 2: This boundary is the interface between the anode 
diffusion layer and the anode rib collector, which is an imper- 
meable wall. Accordingly, all the fluxes in the y direction are 
zero at this boundary: 


F U ômem 3CM 0 dCmy 0 dCwv 0 
With respect to the water transport, the flux due to diffusion dy dy dy 
mode can be ignored if the membrane is in equilibrium with əpı 0 ds 0 15 
liquid water on both sides. Hence, the flux of water crossover ay ay (15) 
Table 1 
Constitutive relations 
Parameters Expressions 
1.41701 — s) — 2.120(1 — s)? + 1.26301 — s), 0 < 8e < 90° 


Capillary pressure Pc = Pg — Pi = 7 COs 6(€/K)° J(s), J(s) = { 


Relative kas? (liquid), krg =(1 — sp (gas) 


permeabilities 


1.4178 — 2.1205? + 1.26357, 


90° < @ < 180° 


Dyiel*s!>, ADL 
E+ EN 
Effective diffusion D = D; ce!(1 — 8)! : Oo, WV, MV), D = , ACL 
coefficients of - ma i i ? A e/(Dyjie's!*) T en/(Dm Nenn) 
species [27] DM,NE ”, MEM 
{ —My,0Rw = Mu Ru. ADL -Mm oŘv, CDL 
General generation Mia Ja a , Mle = je Ip > 
rate of mass in (Mnzo 4 MM) MmoRw — Mu Rm, ACL My,0 OF. 6F8ca MyoRw, CCL 
liquid phase 
My,0 Rw + Mu Ru ADL Mu,0Rw CDL 
General generation Mga = ja ” ? , Mg c = j bp 
rate of mass in gas Mco, 6F + My,oRw + MmRm ACL Mo, AF | Mco, 6Fboa + My,oRw CCL 


phase 


Mole generation rate 
of species 


—Ru ADL | Ru ADL , Rw ADL 
Ru = { oj y Ruv,a = { ee Rwa= 5” 
6F ™ 


Rw ACL 
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e Boundary 4: Similar to Boundary 2, Boundary 4 represents the 
cathode channel rib, at which all the fluxes in the y direction 
are zero: 


3Co, 


dy 


ac 
WV =0 
dy 


e Boundary 5: This boundary represents the inlet of oxygen 
supply and the outlet of water removal at the cathode, at which 
the following boundary conditions are specified: 


an ae 


>, ~=0 16 
Jy By (16) 


’ ’ 


Co, = C, CwH Gy Pe=pp, s=0 (17) 
e Boundaries 3 and 6: These two boundaries are symmetrical, 
respectively, with respect to the middle point of the chan- 
nel width and the middle point of the rib width. Hence, the 


gradients of all the variables in x direction are zero: 


OCM _ OCMV _ dCo, _ OCwv _ 
ax a ax 9 əx 0, əx o o, 
al Ope ðs 
— = 2 = — = 1 
ox 0, Ox 0, Ox ue 


The conditions at the interfaces between those functional lay- 
ers in the ‘sandwiched’ cell are given based on the principle 
that the continuity and mass/species flux balance are required 
at each interface to satisfy the general mass and species conser- 
vation of the entire cell. More details can be found elsewhere 
[27]. 


2.3. Electrochemical kinetics 


On the DMFC anode, the Tafel-like expression is used to 
model the methanol oxidation reaction (MOR) kinetics: 


Y 

ja = Ava (=) exp (Sr ms) (19) 
where the reaction order y is related to the methanol concen- 
tration and is assumed to be zero-order when the methanol 
concentration is higher than a reference value. Otherwise, the 
first-order reaction is specified. 

With respect to the electro-reduction reaction of oxygen on 
the cathode, modified first-order Tafel-like kinetics based on a 
liquid-film covered agglomerate model [27] is used, which gives 


ac F 
se ( ne) ff 2 


: Co,/ku,0 
: ref „O2 
J = Ay,cJ0.0, a." 
02 


R 


where Co, represents oxygen concentration in the gas pores. 
Three factors ky,9,, £1 and z denote the Henry factor to capture 
the effect of dissolving process, the factor in view of the transfer 
resistance in liquid film and the factor in view of the transport 
resistance in the agglomerate, respectively. 


2.4. Current balance and cell voltage 


The protons and electrons produced by the MOR at the anode 
transfer to the cathode through the membrane and the external 


circuit, respectively. The cell current density can be calculated 
by 


SS act Ja dx dy 


(we + wr)/2 os 


Ice = 


To account for methanol crossover, a ‘parasitic’ current den- 
sity, Ip, is used to express the rate of methanol crossover: 


ma 6F {ee NM dx 
j (we + wr)/2 


(22) 


where the molar flux of methanol crossover, Ny, is given by Eq. 
(12): 

At the cathode, it is assumed that both the cell current and the 
‘parasitic’ current are entirely consumed by the oxygen reduc- 
tion reaction (ORR), i.e.: 


SS ccrJe dx dy 


(We + wr)/2 eo 


Icen + Ip = 

For a given anode overpotential, ņa, the cell current den- 
sity, Icey, and the ‘parasitic’ current density can be determined 
from Eqs. (21) and (22), respectively. Then, the cathode mixed 
overpotential, nc, with the effect of methanol crossover, can be 
obtained with the help of Eqs. (20) and (23). Finally, the cell 
voltage can be determined from: 


ô 
Veen = Vo — Na — Ne — Icen (Ress F nen ) (24) 


where Vo, Reontact and « denote the thermodynamic-equilibrium 
voltage of a DMFC, the contact resistance and the proton con- 
ductivity of the membrane, respectively. 


Table 2 

Cell geometric dimensions and operating parameters 

Parameters Symbols Values Units 

Anode diffusion layer thickness Sadl 2.6 x 1074 m 

Anode catalyst layer thickness Sacl 0.2 x 1074 m 

Membrane thickness (Nafion Ô meri 1.3 x 1074 m 
115) 

Cathode diffusion layer thickness Scdl 2.6 x 1074 m 

Cathode catalyst layer thickness bcel 0.2 x 1074 m 

Channel width We 1.0 x 10-3 m 

Rib width Wr 1.0 x 107° m 

Operation temperature T 348.15 K 

Anode inlet pressure p? 1.013 x 1075 Pa 

Cathode inlet pressure Pg 1.013 x 1075 Pa 

Inlet methanol concentration at C in 1000 mol m~? 
anode 

Inlet methanol vapour Chv (eng mol m~? 
concentration 

Inlet oxygen concentration at C Op 7.35 mol m~? 
cathode 

Inlet water vapour concentration Cwv C a mol m~? 
at anode 

Inlet water vapour concentration Cwv 0 mol m~’ 
at cathode 

Inlet liquid saturation at anode sm 1 

Inlet liquid saturation at cathode sin 0 
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3. Results and discussion 


The above-described governing equations for the cell 
geometric dimensions and operating parameters listed in 
Table 2 subjected to the boundary conditions, along with 
the constitutive relations shown in Table 1 and electrochem- 
ical properties listed in Table 3, are solved numerically 
using a self-written code, which was developed based 
on the SIMPLE algorithm with the finite-volume-method 


[27]. 


3.1. Model comparison 


In this section, the cell performance predicated by the present 
model is compared with those predicated by other models. The 
predicted polarization curves at the same operating conditions 
are compared in Fig. 2, where Curve A corresponds to the present 


model, in which transport of methanol in both the liquid and 


Cell voltage (V) 


141 
0.75 
A Two phase model with the effect of non 
equilibriumevaporation/condensation 
B Two phase model without considering 
methanol vapor in gas phase 
C Two phase model with the assumption 
of thermodynamic equilibrium between phases 
0.50 Single liquid-phase model 
0.25 
0.00 
0 1000 2000 3000 4000 5000 6000 


Cell current density (A/m’) 


Fig. 2. Polarization curves predicted by different models with 1.0-M methanol 


operation. 
Table 3 
Physicochemical properties 
Parameters Symbols Values References 
Porosity, permeability 
ADL €adls Kadı (m?) 0.7, 2.0 x 107!? [17] 
ACL €acts Kaci (m?) 0.3, 1.0 x 10-14 
MEM Emem; Kmem (m?) 0.3, 2.0 x 10718 
CCL Eccl; Keel (m?) 0.3, 1.0 x 10-14 
CDL Ecdls Kea (m?) 0.7, 2.0 x 107!? [17] 
Nafion volume fraction 
ACL EN,acl 0.4 
CCL EN cel 0.4 
Diffusivities 
MeOH in water Dy (m? s7!) 6.69 x 107° 23 
MeOH in Nafion Dyn (m? s7!) 4.9 x 107 !0e!2436(0/333)-0/)] 32 
Methanol vapour Dmg (m? s7!) —6.954 x 1076 + 4.5986 x 1078T + 9.4979 x 107!!7° 22 
O% in gas Do, g (m? s7!) 1.775 x 1075(T/273.15)!823 22 
Op in water Doo, (m? s7!) 3.032 x 107° 33 
O% in Nafion Do, N (m? s7!) 1.844 x 10710 33 
Water vapour Dwvg (m? s~!) 2.56 x 1075 (T/307.15)2334 22 
Viscosity of gas phase Mg (kg m7! s7!) 2.03 x 1075 36 
Viscosity of liquid phase wi (kgm7! s7!) 4.05 x 1074 38 
Electro-osmotic drag coefficients of water and methanol nd, H,O 2.5 34 
naM Nd,H»0XM 34 
Evaporation rate constant for water ke (atm! s7!) 5.0 x 1073 21 
Condensation rate constant for water ke (molatm7! s7! cm~?) 5.0 x 1075 21 
Interfacial transfer rate constant for methanol hg (m? s7!) 0.001 24 
Specific interfacial area between liquid and gas Aig (m7!) 105 24 
Proton conductivity in membrane k (Q7! m!) 7.3 el!268(01/298)—(1/7))] 32 
Henry law constant for oxygen kH,o, el(=666/7)+14.D/(RT) 33 
Henry law constant for methanol kum (atm) 0.096 e9.0451 17-273) 22 
Saturation pressure of water vapour log jg pwy (atm) —2.1794 + 0.02953(T — 273 — 9.1837 x 10->(T — 273)? 35 
+ 1.4454 x 1077 (T — 273) 
Saturation pressure of methanol vapour Piy (atm) kH,MXM. 22 
Thermodynamic voltage Vo (V) 1.21 22 
Transfer coefficient of anode Oa 0.5 2! 
Transfer coefficient of cathode Ac 1.0 37 
Anode exchange current density Ay,a J (A m~?) 1.0 x 10° 27 
Cathode exchange current density Ay,c io, (A m~?) 1.14 x 10° 27 
Anode reference concentration cut (mol m~3) 100 22 
Cathode reference concentration cs: (mol m~?) 0.52 27 
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vapour phases is considered and the effect of non-equilibrium 
evaporation and condensation is taken into account to deal with 
the transfer of interfacial species between the phases; Curve 
B corresponds to the model presented in Ref. [27], in which 
the mass transport of liquid methanol and gas carbon diox- 
ide (CO2) is considered, but methanol vapour in the gas phase 
is ignored; Curve C corresponds to the case predicted by the 
present model, but in which the local thermodynamic equilib- 
rium between the phases is assumed. Note that this assumption 
is widely used in a number of fuel cell models reported ear- 
lier [22,23,26]. With the assumption, the gas phase in the anode 
can be treated to be saturated with methanol vapour. The details 
of the model for Curve C are presented in Appendix A. Curve 
D corresponds to the single-phase model, in which mass trans- 
port in the anode is modelled as a pure liquid flow, while in 
the cathode a pure gas flow is assumed. A comparison between 
Curves A and B indicates that the limiting current density pre- 
dicted by the present model is slightly higher than that by 
the one without considering the mass transport of methanol 
vapour. Physically, this means that the added mass transport 
of methanol vapour tends to increase the overall concentra- 
tion of methanol in the anode catalyst layer, and leads to a 
higher limiting current density [22]. Hence, ignoring the mass 
transport of methanol vapour will underestimate the cell perfor- 
mance. 

For the same reason that methanol vapour phase is also 
considered in predicting Curve C, the limiting current den- 
sity in Curve C is also higher than that in Curve B. However, 
a rather large difference can be found when a comparison 
is made between Curves A and C despite the fact that the 
transport of methanol vapour is considered in both corre- 
sponding models. This difference is attributed mainly to the 
fact that the model in predicting Curve A takes account 
of the effect of non-equilibrium evaporation and conden- 
sation, whereas the one in predicating Curve C has the 
assumption of a local thermodynamic equilibrium between the 
phases. 

To illustrate which model gives a more realistic predic- 
tion, Curves A and C are now compared with Curved D, 
which is predicted by a single-phase model. As shown in 
Fig. 2, the limiting current density in Curve C is much higher 
than that in Curve D. However, the study has shown that the 
single-phase model over-predicts cell performance, because the 
presence of gas bubbles in the anode porous region and in the 
flow channel increases the resistance of methanol transfer to 
the catalyst layer and thereby results in an increase in mass- 
transport polarization [20,39]. Therefore, the models that have 
the thermodynamic-equilibrium assumption may overestimate 
the two-phase mass transport to give an unreasonable predic- 
tion of cell performance. By contrast, it is seen from Fig. 2 that 
the limiting current density in Curve A is smaller than that in 
Curve D, meaning that the present model with the effect of non- 
equilibrium evaporation and condensation gives a more realistic 
prediction of cell performance. In the following sections, all 
the results to be reported are predicted by the present model 
with the effect of non-equilibrium evaporation and condensa- 
tion. 


3.2. Velocity field and concentration distributions of 
methanol at different phases 


To illustrate the liquid—gas two-phase flow behaviour, the 
velocity fields for both liquid and gas phases are shown in Fig. 3. 
As shown in Fig. 3a, the liquid phase moves from the anode chan- 
nel to the catalyst layer as the result of the water and methanol 
consumption by the MOR and by the water and methanol 
crossover through the membrane to the cathode. By contrast, 
Fig. 3b shows that the gas phase in the anode porous region 
moves in the opposite direction of the liquid phase, as CO2 is 
generated in the catalyst layer and moves towards the flow chan- 
nel region. The liquid—gas two-phase counter-flow behaviour in 
the anode porous region shown in Fig. 3 illustrates the com- 
plexity in simultaneously transporting the reactants (methanol 
and water) towards the catalyst layer and removing the product 
(CO2) from the porous region. 

In order to indicate the two-phase mass-transport behaviour 
in the DMFC anode, the distribution of liquid methanol concen- 
tration, methanol vapour concentration and gas-void fraction 
are depicted in Fig. 4. Clearly, as shown in Fig. 4a, methanol 
concentration is highest at the anode-channel—diffusion-layer 
interface, and then decreases from the region near the rib region 
and in the catalyst layer as the result of mass-transfer resistance 
and mass consumption by the MOR in the catalyst layer. The 
gradient of methanol concentration implicitly indicates the dif- 
fusive transport of methanol, while the motion of liquid in porous 
anode as shown in Fig. 3a implies the convective transport of 
methanol. As a result, the transport of methanol in the liquid 
phase is determined by the superposition of two mechanisms: 
molecular diffusion and macroscopic convection. 

The distribution of the methanol vapour concentration is pre- 
sented in Fig. 4b. The concentration is higher near the channel 
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Fig. 3. Distributions of velocity vector of (a) liquid and (b) gas anode porous 
region (anode overpotential: 0.5 V). 
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Fig. 4. Distributions of (a) methanol concentration in liquid phase through anode electrode and membrane, (b) methanol vapour concentration in gas phase and (c) 


gas-void fraction in anode porous region (anode overpotential: 0.5 V). 


region and becomes lower in the catalyst layer as the result of 
the MOR there. Fig. 4c gives the distribution of the gas-void 
fraction, œ, in the anode porous region. (Note that the value of œ 
is obtained by subtracting liquid saturation from unity.) Clearly, 
the gas coverage level is higher in the region close to the catalyst 
layer than in the region close to the channel due to the generation 
of CO7 as the result of the MOR in the catalyst layer. Moreover, 
the gas tends to accumulate in the impermeable region under the 
rib and this results to a higher level of gas coverage in the region 
under the rib. An increase in gas coverage in the porous medium 
will reduce the space available for methanol solution, which, in 
turn, will lead to an increase in the mass-transfer resistance of 
methanol from the channel to the reaction catalyst layer. 


3.3. Analysis of mass transport of methanol 


Next, we consider the two-phase mass transport of methanol 
in the DMFC anode. The total flux of methanol in the anode 
diffusion layer is comprised of two contributions, namely: 
methanol transport in a liquid solution and methanol vapour 
transport in the gas phase. The variation in current density in 
the catalyst layer along the channel and rib width (x direction 
in Fig. 2) for three different cases is shown in Fig. 5. Case A 
corresponds to the model that only considers methanol diffusive 
transport in a liquid phase; Case B considers methanol diffu- 
sion and convection transport in the liquid phase; and Case C 
corresponds to the model that considers methanol diffusion and 
convection in both the liquid phase and the gas phase. From 
Fig. 5, clearly, the current density for Case C is higher that for 
Case B owing to the contribution of methanol vapour transport 
in the gas phase. It also can be seen that the current density 
for Case B is a slightly higher than the value for Case A due 
to the fact that the convective motion of liquid solution tends 
to decrease the resistance of methanol transport to the catalyst 
layer, as discussed earlier. 


The limiting current density corresponding to the mass- 
transport limitation of methanol as the result of each of the 
three transport mechanisms, A-C, is shown in Fig. 6 for 0.5, 
0.75 and 1M methanol feeding concentrations. Note that in 
Fig. 6 the scale of the limiting current density was set to begin 
with 1000 A m~? to give prominence to the differences between 
Cases A-C. It should be recognized that the difference in the 
limiting current density between Cases C and B corresponds to 
the mass flux of methanol due to transport in the vapour phase, 
designated by I; the difference between Cases B and A corre- 
sponds to the mass flux of methanol due to convection in liquid 
phase, designated by II; and II in Fig. 6 represents the mass 
flux of methanol due to diffusion in liquid phase, i.e., Case A. 
As can be seen from Fig. 6, for a 1.0M methanol feeding con- 
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Fig. 5. Limiting current densities corresponding to mass flux of methanol as 
result of (A) methanol diffusion in liquid solution only; (B) methanol diffusion 
and convection in liquid solution; (C) methanol diffusion and convection in both 
liquid and gas phase. 
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Fig. 6. Limiting current densities corresponding to mass flux of methanol as 
result of different transport mechanisms for different methanol concentrations. 


centration, the contribution by methanol transport in the vapour 
phase and methanol convective transport in the liquid phase, 
respectively, take about 6% of the total limiting current density. 
The rest 88% is contributed by methanol diffusive transport in 
the liquid phase. For all the cases studied, the results indicate 
that the mass flux of methanol is predominated by the methanol 
diffusion in the liquid phase. Methanol convection in the liq- 
uid phase and methanol vapour transport in the gas phase only 
contribute a relatively small portion of the total mass flux. How- 
ever, it should be pointed out that the results presented in Fig. 6 
were obtained for a DMFC with a parallel flow-field, in which 
the convective motion of liquid is relatively weaker. For a ser- 
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pentine flow-field, however, the pressure difference between the 
two adjacent channels will force the liquid convection through 
the diffusion layer and this induces stronger forced-convection 
[17,31]. As aresult, methanol convection in the liquid phase may 
become more important than methanol diffusion in the liquid 
phase. 

Another point that should be mentioned about the results in 
Fig. 6 is that methanol vapour transport in the gas phase is influ- 
enced by the total gas-void fraction in the porous region. In this 
study, the gas-void fraction in the anode porous region, as shown 
in Fig. 3c, is relatively small because of the assumption that 
no gas covers the anode-channel—diffusion-layer interface. This 
assumption can be justified only when gas bubbles in the anode 
channel are removed quickly by the high-speed liquid flow. Oth- 
erwise, coverage of the interface with gas bubbles will result in a 
much higher gas-void fraction inside the diffusion layer and the 
catalyst layer. As a result, the contribution of methanol vapour 
transport in the gas phase will play an important role. The influ- 
ence of the interfacial gas coverage will be discussed in the 
following section. 


3.4. Effect of gas coverage level at the 
anode-channel-diffusion-layer interface 


In the DMFC anode, it is understood that the presence of 
gas bubbles accumulated in the flow channel will reduce the 
contact surface area between the liquid methanol and diffusion 
layer, thus blocking methanol transfer from the channel to the 
diffusion layer. Fig. 7 shows the distribution of the gas-void 
fraction in the porous anode when the interfacial gas coverage 
level is at æ =0.2, 0.4 and 0.6. Generally, the gas-void fraction 
is higher in the region close to the catalyst layer and in the 


0.001 


0.00028 0.00028 
Y Y 


(b) (c) 


Fig. 7. Distributions of gas-void fraction in anode porous region for different interfacial gas coverage levels: (a) @arc/ADL = 9.2; (b) @AFCc/ADL =0.4; and (c) 


Q@AFC/ADL = 0.6 (methanol concentration: 1 M; anode overpotential: 0.5 V). 
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Fig. 8. Distributions of local current density along channel and rib width for 
different interfacial gas coverage levels (methanol concentration: 1.0 M; anode 
overpotential: 0.5 V). 


region under the rib than in the region close to the flow channel. 
It can be seen that the overall gas-void fraction increases with 
increasing gas coverage level at the interface. This result explains 
why accumulation of gas at the AFC—ADL interface hinders the 
product gas CO2 moving out of the diffusion layer and results 
in accumulation of gas bubbles in the porous region. 

The variation in the local current density along the channel 
width and the rib width for different interfacial gas coverage 
levels is presented in Fig. 8. It is clear that the overall local cur- 
rent density is lowered with an increasing level of interfacial gas 
coverage. The effect of the degree of interfacial gas coverage 
on cell performance is presented in Fig. 9. As can be seen, the 
limiting current density decreases with increasing the interfacial 
gas coverage. It is worth mentioning that, in practice, the level 
of interfacial gas coverage cannot be arbitrarily specified, as it 
is influenced by many factors such as the anode flow rate, the 
geometry of flow channel, and the interfacial properties of the 
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Fig. 9. Effect of interfacial gas coverage on cell performance with 1.0M 
methanol operation. 


gas-diffusion layer. Moreover, the extent of interfacial gas cover- 
age also varies along the flow direction due to the accumulation 
of gas along the channel. 


3.5. Effect of cell dimensions 


The mass transport inside the DMFC is influenced by many 
parameters of the structural details, such as the channel width 
(wch) and the channel-to-rib width ratio (C/R). In order to 
improve cell performance, it is essential to understand the para- 
metric effects on two-phase mass transport and cell performance. 


3.5.1. Effect of channel width 

To study the effect of channel width, the C/R ratio was fixed at 
unity. It is should be recognized that reducing the channel width 
while fixing the C/R ratio is equivalent to having more channels 
in the parallel flow-field. The distribution of the gas-void frac- 
tion and the methanol concentration at the ADL—ACL interface 
along the channel and the rib width (in the dimensionless form, 
X = x/wWen) for different channel widths is presented in Fig. 10. 
As shown in Fig. 10a, the overall level of the gas-void fraction 
at the ADL—ACL interface decreases and becomes a more uni- 
form distribution as the channel and rib widths decrease. It is 
understood that the produced CO2 in the region under the rib 
has to be transported across the thickness of the ADL and the 
width of the rib to the channel. When the channel and the rib 
widths decrease while the thickness of the diffusion layer does 
not change, the shorter width of rib assists removal of the gas in 
the region under the rib. This explains why significant changes 
in the gas coverage level occur in the region under the rib as the 
rib width changes. By contrast, only small changes are observed 
in the regions under the channel. The decrease in gas coverage 
results in a larger area for the transport of methanol solution. 
It is seen in Fig. 10b that the overall methanol concentration is 
more uniform and higher for a smaller rib width as the result 
of the short transfer distance from the channel to the region 
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Fig. 10. Distributions of (a) gas-void fraction and (b) methanol concentration at 
ADL -ACL interface for various channel widths (methanol concentration: 1.0 M; 
anode overpotential: 0.5 V). 
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under the rib and decreased gas blockage. The uniform distribu- 
tion of methanol in the catalyst layer is favourable for sufficient 
utilization of the MEA area under the rib region. 

Reducing the channel and the rib width leads to higher 
methanol concentration and lower gas coverage in the region 
under the rib, resulting in higher and more uniform current den- 
sity distribution. As a result, the limiting current density of the 
cell increases, as can be seen in Fig. 11. It must be pointed out 
that although a smaller channel is preferred if the C/R ratio is 
constant, one has to consider the fact that a smaller channel width 
may result in a higher flow resistance in the channel. Moreover, 
a smaller channel may also result in transient gas bubble hold- 
up in the flow channels, thereby severely blocking the methanol 


supply [20]. 


3.5.2. Effect of channel-to-rib width ratio 

The effect of the C/R ratio is studied by keeping the total width 
of the channel and the rib constant. The distributions of local 
gas-void fraction and methanol concentration at the ADL-ACL 
interface for different C/R ratios are shown in Fig. 12. As the 
C/R ratio increases, the overall gas-void fraction decreases and 
becomes more uniform due to the decrease in transport distance 
for gas removal under the rib, which thus offers more transfer 
space for methanol solution. Accordingly, the methanol concen- 
tration increases and becomes more uniform in distribution with 
increasing C/R ratio, as shown in Fig. 12b. The results indicate 
that a larger C/R ratio leaves a larger accessible catalyst area 
and also improves the removal of products, resulting in a more 
uniform distribution of reactants and a more effective utilization 
of the MEA surface. 

The effect of the C/R ratio on cell performance is presented 
in Fig. 13. The results show that the limiting current density 
increases with increasing C/R ratio, due mainly to the increased 
transfer area exposed to the open channel. It should be rec- 
ognized, however, that an increasing C/R ratio may lead to a 
decrease in the contact area between the rib and electrode and the 
electron collector which, in turn, may cause an increase in cell 
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Fig. 11. Effect of channel width on cell performance (methanol concentration: 
1.0M). 
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Fig. 12. Distributions of (a) gas-void fraction and (b) methanol concentration 
at ADL-ACL interface for different C/R ratio (anode overpotential: 0.5 V). 
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Fig. 13. Effect of channel-to-rib width ratio on cell performance (methanol 
concentration: 1.0 M). 


internal resistance. A larger C/R ratio also means that electrons 
have to transfer over a longer distance at the in-plane direction in 
the diffusion layer to the current-collector, resulting in a larger 
Ohmic loss. This aspect must be taken into consideration in the 
optimum design of the DMFC. 


4. Conclusions 


In this work, a two-phase mass-transport model for a 
liquid-feed DMFC has been developed. It considers the trans- 
port of methanol/methanol vapour, and water/water vapour at 
the anode, and the transport of oxygen at the cathode. In 
particular, the model takes into account the effect of non- 
equilibrium evaporation and condensation, which eliminates the 
thermodynamic-equilibrium assumption employed in other pre- 
vious two-phase mass-transport models. It has been shown that 
the thermodynamic-equilibrium assumption between liquid and 
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vapour phases will lead to an unreasonable prediction of cell 
performance. The parametric study using the model developed 
in this work indicates that: (1) the gas coverage at the anode- 
channel—diffusion-layer interface has a significant influence on 
the gas-void fraction inside the anode porous region, thereby 
increasing the mass-transfer resistance of methanol and lower- 
ing cell performance; (2) the geometric dimensions of the cell 
structure, including the rib width and channel-to-rib width ratio, 
greatly influence the uniformity of reactant distributions, the 
current density profile, and the overall cell performance. 
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Appendix A. The model with the assumption of 
thermodynamic equilibrium 


The assumption that a local thermodynamic equilibrium pre- 
vails between the liquid and gas phases implies that 


CMV. = = (al) 


where the methanol vapour pressure can be obtained by Henry’s 
law: 


sat M H20 


Puv = KH,MXM,1=kH,M 


CM. (a2) 
Pl 


The mass conservation of methanol can now be obtained by 
superposition of methanol transport in the liquid and gas phases 
to give 


y. [(- Kky 
Hı 


Kkrg eff 
+V- ||- 7 V Pg,a | CMV, g,a — DMv,g V CMV.g,a 
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va) CM. — DEV Oma 


= Rua + ŘmMv.g,a (a3) 


Combining Eqs. (al)—(a3), we obtain 


Kku Kkrg 
V. = Vpiat | — WY Pg,a Cla 
Hı Mg 


— (Dif + Div eV Cua} = Ruor (a4) 


where the term Ryor represents the molar consumption rate 
by the MOR in the anode catalyst layer only. Solving Eq. (a4), 
along with Eqs. (3)—(10), yields the results for the case in which 


it is assumed there is thermodynamic equilibrium between the 
phases. 
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